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Abstract: A recent experimental study [1] has demonstrated the alignment of phase 
separated domains across hundreds of bilayer units in multicomponent stacked lipid 
bilayers. The origin of this alignment is the interlamellar coupling of laterally phase 
separated domains. Here, we develop a theoretical model that presents the energetics 
description of this phenomenon based on the minimization of the free energy of the system. 
Specifically, we use solution theory to estimate the competition between energy and 
entropy in different stacking configurations. The model furnishes an elemental phase 
diagram, which maps the domain distributions in terms of the strength of the intra- and 
inter-layer interactions and estimates the value of inter-layer coupling for complete 
alignment of domains in the stacks of five and ten bilayers. The area fraction occupied by 
co-existing phases was calculated for the system of the minimum free energy, which 
showed a good agreement with experimental observations. 
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1. Introduction 

Phase separation and domain formation in multicomponent lipid bilayers containing immiscible 
liquid phases is a well-known phenomenon. In this regard, ternary lipid mixtures consisting of 
cholesterol (Choi), sphingomyelin (SM), and unsaturated phospholipids are of considerable interest as 
models for lipid rafts, which are thought to be relevant in a variety of cell-surface signaling in 
biological membranes [2-5]. 

These lipid mixtures, below their miscibility transition temperatures, readily phase-separate into 
two co-existing liquid phases: (1) liquid-ordered (L 0 ) and (2) liquid-disordered (Ld) [6-9]. L 0 domains, 
or the so called lipid rafts, predominantly consist of Choi and SM. Ld domains are enriched in 
unsaturated phospholipid. Recently, we investigated these raft-forming ternary mixtures as stacked 
lipid bilayers (lipid multilayers), which led to the interesting observations of the long-range interlayer 
alignment of laterally phase separated L 0 domains [1]. We found that the domains formed in each 
bilayer become aligned across several lamellae in a stacked lipid bilayer. 

The demonstrated bio-epitaxial film is of considerable technological interest for biomimicry. 
For instance, thylakoid membranes of photosynthetic cyanobacteria [10] or plant chloroplasts [11], and 
electrocyte cells in electric eels [12] contain multilamellar membranous structures and their mimicry 
can be exploited in energy harvesting technologies. Potential applications of such membranes in 
photonics and sensing are also significant [13-15]. 

In the experimental investigation presented in [1] various compositions were examined including 
equi-molar SM and DOPC mixtures with different concentrations of cholesterol in the range of 
10-40 mol% of the total solution [Cholx/SMl/DOPCl; x = 10%-40%]. The domain formation and 
their alignments were observed in all of such compositions. Their alignments were examined using a 
combination of X-ray diffraction, fluorescence confocal microscopy, and atomic force microscopy 
measurements. In this study, we also carried out additional experiments to validate the theoretical 
model, which is the main subject of the current paper. 

Examples of the formation of the aligned domains in lipid multilayer and their epifluorescence and 
atomic microscopy (AFM) images are shown in Figure la,b. As can be seen in Figure lcl-c3, 
entropically a multilayer system prefers a random distribution of the domains. Alignment of the 
domains requires the existence of a strong coupling between bilayers which can overcome the entropy. 
It was suggested that the origin of this coupling is the energy penalty associated with the tension 
between the different water networks at the vicinity of the different domains [1]. It was suggested that 
the origin of this coupling is the energy penalty associated with the tension between the different water 
networks at the vicinity of the different domains [1]. 

To analyze this coupling and estimate its value, here we investigated the competition of the energy 
and the entropy in multicomponent lipid multilayer systems with varying number of lamellae. 
We derived an elemental phase diagram that can explain the relation between the inter- and intra-layer 
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coupling parameters. The phase diagram can predict the production of the multilayers with different 
domain distributions (different phases). It can further estimate the intra-layer coupling value for the 
complete alignment of the multilayers. The results are presented for the two cases of 5 and 10 number 
of lamellae. The model calculation can also predict the areal fraction of the multilayer surface adopted 
by different phases. Presenting an example of this calculation, we show that the theoretical result is in 
good agreement with the corresponding experimental data. 

Figure 1. (a) Epiflorescence image of the top surface of a multi component lipid 
multilayers doped with 0.03 mol% of Rhodamine which partition preferebly with Ld 
domains. Scale bar is 100 |um; (b) Atomic force microscopy image of the multilayer cleraly 
shows the aligned domains; (cl,c2) Entropy favors the randome distribusion of the 
domains in different bilayers of the multilayer as shown in panel (cl); However, the 
inter-layer couling competes with entropy to align the domains. More coupling strengh 
lead to more complete alignmnet (c2,c3). 




2. Results and Discussions 

Although the approach to the equilibrium for the domain dynamics in multilamellar lipid 
membrane, specifically at low humidity, is slow, membrane multilamellae ultimately do reach to a 
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point that no further changes are observed at experimental time scales (several days) [1]. Under these 
conditions of quasi-equilibration, applying energetic considerations becomes instructive, such as the 
model described below. 

To model the experimental system of multicomponent lipid mulilayer, we begin by considering 
each bilayer as a consolidated unit. In the other words, we assume that the domains, formed via 
intralamellar phase separation, are pre-aligned across the two leaflets of each bimolecular lamellar 
sheet (monolayers). Exactly similar to the experimental report in reference [1], we consider the 
composition of Choi, SM and l,2-dioleoyl-sn-glycero-3-phosphocholine (DOPC) as the unsaturated 
lipid. Next, considering nearest-neighbor intra-layer interaction between the species, which is positive 
(unfavorable) between DOPC and Choi and negative (favorable) between SM and Choi [2,16], we 
assume SM and Choi represents one component and model our system as a pseudo-binary mixture [17] 
with component 1: DOPC and component 2: SM/Chol. These two components undergo phase 
separation to make two phases DOPC-rich and L 0 : SM/Chol-rich. To account for the small 
non-ideal compositional variation in each layer, a normal distribution with mean value equal to the 
molar ratio of the composition and a small standard deviation of 0.02 is assumed. 

The total free energy, F 9 of the multilamellar membrane can be approximated using the mean-field 
description of the regular solution theory. Within this framework, (1) assuming the random mixing 
approximation of Bragg- Williams type for each layer [18,19], (2) introducing a coupling term for the 
interaction between bilayers, (3) adding the additional entropy term due to the ordering multiplicity of 
the multilayer, and finally (4) summing over all lamellae, we can write: 

F = X 'j{ e jSj) 

0 1 ^ 1 +0 2 ^ 2 j=Phase number ' ^ 



Where 



TV 

e j = S [mj (i - <Pn ) + \j - <p MJ ) 2 ) (2) 



(3) 



Where "F 9 is the total free energy of the multilayer in units of K B T (thermal energy), and "f is the mol 
fraction of each "phase". Here, the term "phase" refers to a multilamellar region with specific local 
composition in each layer (e.g., Figure 2a) (Note that this definition of "phase" differs from L 0 or 
La phase). 
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Figure 2. Simulation results for a 10-layer sample with X — 10%. (a) Three phases are 
predicted for % = 2.3 and A = 0.7: phase 1: aligned Ld, phase 2: unaligned, and phase 3: 
aligned L 0 . The local molar ratio of DOPC (cp) for each layer in each phase (q>#) is shown in 
its corresponding box. The brighter the area, the higher the concentration of DOPC. 
The width of each box is proportional to the surface area of each domain. Each phase 
occupies "t mole%", and "A area%" of the multilayer top surface; (b) Image analysis of 
several experiments in multilayer atX= 10% is summarized in this plot (variation bars are 
shown). The area occupied by the bright, gray, and dark domains (representing the 
fluorescent intensity in the experiment) is in good agreement with the area predicted for 
phases 1, 2 and 3 in plot a, respectively; (c) An example of a phase separated lipid 
multilayer: Gray regions indicate the unaligned domains. Each specific color of gray 
represent one "phase". When there are various gray regions in the multilayer, the system is 
called multi-phase. When there is no gray region, the system is composed of only two 
phases (dark and bright). In this case we have a "two-phase" multilayer in which the 
domains are completely aligned. Note that this picture is chosen as the gray regions can be 
clearly seen so we can clarify the concept of the "phase" in our model. Usually, multilayers 
do not consist of such wide regions of gray area and very small area fraction is adopted by 
gray region as it is shown in Figure la. 
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As it is described in Figure 2c, for this terminology, experimentally, different phases are 
distinguished by their corresponding intensities observed in the epifluorescence measurements. 
The intensity of each phase, as seen from the top surface, is the superposition of the intensity of all 
layers under that region. For instance, in the case of only two distinguishable intensities, dark and 
bright, the multilayer is represented in terms of two phases. 

"<Di" and "® 2 " are the input values of the model and indicate the average mole value of the first and 
second component in the media, respectively. In other words, they are the mean composition (molar 
ratio) of component 1 (DOPC) and component 2 (SM/Chol) in the primary total stock solution, "a" 
and "02" are the molecular surface area of each component. "A" represents the total top surface area of 
the multilamellar sample. "E" and "S" refer to the energy and entropy parts of the total free energy in 
each phase, respectively. "N" is the number of lamellae. "9,/' is the local composition of one 
component (e.g., DOPC) in phase j and layer /. "%" is the degree of non-ideality expressing the strength 
of the nearest-neighbor interaction in a given layer. For simplicity, we refer to it as intra-layer coupling 
parameter. We assume % is similar for all the layers. "A" is the inter-layer coupling parameter, which is 
assumed to be constant (i.e., Ay = A). Each phase j may consist of several repeating (p^s, which can be 
grouped together (e.g., in the first phase in Figure 2a two <p = 0.80 make a group). "&" is the number of 
such groups in that phase (e.g., in the first phase in Figure 2a k = 2). "n(£ 9 j)" is the number of 

members in each group (e.g., in Figure 2a, «(1,1) = 2 and «(2,1) = 2). In Equation 2, the first term 
corresponds to the nearest-neighbor interaction in each layer. The second term is the energy penalty 
due to the local compositional mismatch between the adjacent bilayers. In Equation 3, the first two 
terms are the solution entropy associated with each layer. The last term is due to the entropy increase 
by the additional multiplicity due to repeating cp^s in each phase j. This term increases significantly 
with the number of lamellae, aj is taken 0.69 nm for all compositions [20]. ci2 is different for different 
concentrations of Cholesterol. ForX= 10%, 15%, 20%, 30% and 40%, values are taken as 0.48, 
0.43, 0.40, 0.40, 0.39 and 0.39 nm 2 , respectively [21]. To find the most stable 9^ and the 
corresponding tjs 9 F is minimized with respect to all variables under the following constraints: 

j J 

Where "® " is the average composition of one component (e.g., DOPC) in layer /. 

To further clarify the concepts of "phase" and "q>,/' let us consider a multilayer which is observed 
from the top surface such as the one shown in Figure 2c. In fluorescence imaging, the intensity in each 
point of the top surface of the multilayer is the superposition of intensities of all layers on that point. 
As a result, since the fluorescent molecules partition with only one type of lipid, the different colors 
(grays with different intensities) that can be seen on top surface of the multilayer, represent different 
columnar series of local lipid composition of all layers, which we call them different "phases". 
For example, if we could visually distinguish three distinct colors on top of the total area of one 
multilayer, we would have three phases of that multilayer. Using the above equations, we can also 
calculate the number of phases in each multilayer. For example, Figure 2a which shows the result of 
the calculation performed for a multilayer consisting of ten lamellae (N = 10), contains three phases. 
A phase difference less than the variance of compositional fluctuation is ignored in the model. 
Each phase is shown by a column with 10 layers (10 boxes). The composition of DOPC in each phase 
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and each layer is calculated and written in the corresponding box. This local composition of one 
specific lipid component in one phase is defined as <p# (/ = layer number and j = phase number). 
Other than the values of all q># in the multilayer system, the "t mole %", and consequently the "area 
fraction %" of the multilayer top surface for all distinct phases, is also included in the outcomes of the 
model. These values give us the chance to compare the results of our model with the experiment as we 
can experimentally measure the area fraction of distinct phases in different multilayers. 

To compare and verify the theoretical results with the experiment, several samples with similar 
composition were prepared following the method presented in [1]. In brief, the stack of the lipid 
mixture was incubated for two days at 50-60 °C followed by decreasing the temperature to room 
condition (-22 °C). Imaging and Epifluorescence measurements were performed using an inverted 
epifluorescence microscope (Nikon TE200E, Technical Instruments, San Francisco, CA, USA) 
equipped with a high spatial resolution CCD camera. The image analyses for estimating the surface 
area of the bright, dark and gray regions were performed using NIH Image J. Figure 2b shows a good 
agreement between the theoretical predictions and the experimental data. 

Increasing A decreases the number of phases and the area of the unaligned regions. A t h res is the 
lowest coupling strength that gives rise to the complete alignment (two fully aligned phases). In the 
phase diagram presented in Figure 3a, this value is calculated for two cases of TV = 5 and N = 10 at 
different % values. The value of % depends on the composition of the system. For larger % values, larger 
coupling strength is required for complete alignment. Also, below a critical value for x, %cnt 9 the 
multilayer does not phase separate and remains completely mixed, jfcn* is close to 2 and increases 
slightly with the coupling parameter. Two curves jfcn* and A t h res divide the phase diagram into three 
regions of mix, multi-phase, and two-phase. For comparison Athres for N = 5 is also plotted (Figure 3a). 
Figure 3b shows that Athres increases with the degree of lamellarity. This is reasonable since the 
membrane multilayers with greater number of lamellae would require stronger interbilayer coupling 
strength to overcome the entropy. 

In brief, Figure 3, clearly indicates that by increasing the number of lamellae, stronger coupling is 
required to align the domains, as expected. This is indeed in agreement with the experimental results 
presented in [1] that a multilayer with the number of lamellae more than -1000 cannot be aligned 
under normal conditions. 

Extending the model to a larger number of lamellae is computationally too expensive. Nevertheless, 
we expect similar phase diagram for larger number of lamellae once the threshold coupling parameter 
is accurately estimated. 

It is important to mention that, prior to the recent experimental data reported in [1], the alignment of 
the domains in lipid multilayers was assumed to be impossible as entropy cannot trivially let formation 
of such ordering in the distribution of the domains. Moreover, previous reports of the alignment of the 
domains in multilayers with low number of lamellae (<5) attached to the substrates showed that the 
bilayers do not couple strongly to align the domains [22-24]. We believe that the reason of the 
unsuccessful alignment of the domains in the reported experiment is the influence of the substrate, 
which is most pronounced for the first few layers close to the substrate. Thus, we may conclude that 
although the coupling is strong, it cannot overcome the interaction between the substrate and the first 
few bilayers. Such an interaction is not taken into account in our model and the layers are assumed to 
be away from the substrate by at least a few layers. 
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The simulation was written and run using Matlab R2010a and Image analyses were performed 
using NIH Image J. 

Figure 3. The phase diagram defining different regions of "mix", "two phase", and 
"multi phase" versus in-plane % and cross-plane A coupling parameters. At and above A t hres 
line, there are only two phases, "aligned Ld" and "aligned L 0 ". Under A t h re s line, the system 
is multi-phase. Under j^ rih the system does not undergo phase-separation and will be fully 
mixed, c, An example of dependency of Athres to the number of lamellae for % = 2.5. 
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3. Conclusions 

We presented a theoretical model which describes the phase separation in multicomponent lipid 
multilayers and the coupling involved in the alignment of the domains across several lamellae. 

This work elucidates the phenomenon recently reported regarding the long-range interlayer 
alignment of the intra-layer domains in lipid multilayers [1]. The stacked lipid bilayers with such 
aligned domains have fundamental biological importance [10,11,25] with potentials for practical 
applications [13-15,26] which necessitate the detailed theoretical analysis of the subject. 

We introduced an elemental phase diagram for multilayers with up to 10 bilayers. The value of the 
interlayer coupling for complete alignment was further calculated for the membranes with different 
intra-layer coupling parameters. 

The surface area fraction of the multilayers occupied by different phases in different situations can 
be calculated with this model. To verify the model, experimental data were generated for a case study 
and the results were compared with those of the model calculations, which showed a good agreement. 
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